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Particle  simulation  methods  are  used  to  address  issues  related  to  recent  results  of  the  NRL  laser  experiment.  The 
simulations,  which  include  all  the  electron  and  ion  scale  lengths,  examine  dense  plasma  clouds  expanding  across  a 
magnetic  field  into  a  vacuum  under  various  conditions.  The  wavenumbers  of  the  observed  short  wavelength  instability, 
identified  as  the  lower  hybrid  drift  instability,  are  found  to  be  independent  of  the  magnetic  field  (in  agreement  with 
experiment)  and  the  ion  to  electron  mass  ratio,  consistent  with  a  simple  model  for  the  expansion.  In  addition,  differences 
in  the  instability  due  to  the  initial  velocity  distribution  and  the  angular  spread  of  the  expansion  are  studied.  Nonlinear 
phenomena,  such  as  possible  flute  tip  bifurcation,  bending  and  separation,  have  also  been  observed  in  the  calculations 
and  are  related  to  similar  effects  seen  in  the  experiments.  However,  detailed  analysis  of  the  simulations  show  that 
there  is  no  short  wavelength  cutoff  to  the  instability,  contrary  to  the  predictions  of  linear  theory,  indicating  that  the 
wavelength  scaling  is  spurious.  Possible  physical  and  numerical  reasons  for  this  behavior  are  examined,  and  it  is  shown 
that  another  instability,  the  electron  cyclotron  drift  instability  can  still  be  excited  at  shorter  wavelengths.  The  overall 
consequences  of  the  results  on  the  suitability  of  these  methods  for  modeling  the  laser  experiment  are  discussed. 
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Application  of  Particle  Simulations  to  the  NRL  Laser  Experiment 

Dan  Winske 

Los  Alamos  National  Laboratory 

Abstract 

Particle  simulation  methods  are  used  to  address  issues  related  to  recent  results  of 
the  NRL  laser  experiment.  The  simulations,  which  include  all  the  electron  and  ion  scale 
lengths,  examine  dense  plasma  clouds  expanding  across  a  magnetic  field  into  a  vacuum 
under  various  conditions.  The  wavenumbers  of  the  observed  short  wavelength  instability, 
identified  as  the  lower  hybrid  drift  instability,  are  found  to  be  independent  of  the  magnetic 
field  (in  agreement  with  experiment)  and  the  ion  to  electron  mass  ratio,  consistent  with 
a  simple  model  for  the  expansion.  In  addition,  differences  in  the  instability  due  to  the 
initial  velocity  distribution  and  the  angular  spread  of  the  expansion  are  studied.  Nonlin¬ 
ear  phenomena,  such  as  possible  flute  tip  bifurcation,  bending  and  separation,  have  also 
been  observed  in  the  calculations  and  axe  related  to  similar  effects  seen  in  the  experiments. 
However,  detailed  analysis  of  the  simulations  show  that  there  is  no  short  wavelength  cutoff 
to  the  instability,  contrary  to  the  predictions  of  linear  theory,  indicating  that  the  wave¬ 
length  scaling  is  spurious.  Possible  physical  and  numerical  reasons  for  this  behavior  are 
examined,  and  it  is  shown  that  another  instability,  the  electron  cyclotron  drift  instability 
can  still  be  excited  at  shorter  wavelengths.  The  overall  consequences  of  the  results  on  the 
suitability  of  these  methods  for  modeling  the  laser  experiment  are  discussed. 


1.  Introduction 


Plasmas  expanding  (or  being  compressed)  at  subAlfvenic  speeds  perpendicular  to  a 
magnetic  field  develop  flute  instabilities  at  the  plasma-magnetic  field  interface.  The  effect 
has  been  known  for  many  years  in  theta-pinch  implosions  [Keilhacker  et  al.,  1974]  and  has 
also  been  observed  in  laser  produced  plasmas  [Okada  et  al.,  1981].  More  recently,  the  same 
phenomena  has  been  seen  in  expanding  barium  clouds  in  the  magnetotail  [Bernhardt  et 
al.,  1987]  and  in  laser-plasma  interactions  [Ripin  et  al.,  1987].  In  the  laser  experiments  the 
flute  structures  exhibit  very  interesting  behavior:  the  outer  tips  are  freestreaming  (i.e.,  not 
slowed  down  by  the  magnetic  field),  sometimes  seem  to  bifurcate,  and  at  late  times  often 
have  a  spiral  pattern. 

A  theoretical  explanation  for  the  laser  observations  has  been  put  forward  by  Huba  et 
al.  [1987],  extending  earlier  work  for  the  barium  releases  [Hassam  and  Huba,  1987].  On 
the  basis  of  a  modified  set  of  MHD  equations,  they  suggest  that  the  instability  which  gives 
rise  to  the  flutes  is  an  unmagnetized  ion  Rayleigh- Taylor  interchange  mode  that  arises  due 
to  the  slowing  of  the  plasma  as  it  expands  against  the  magnetic  field.  The  theory  predicts 
a  growth  rate  that  increases  linearly  with  k,  the  wavenumber,  with  L l/2,  Ln  being  the 
density  scale  length,  and  g 1/2,  where  g  is  the  plasma  deceleration.  Alternatively,  Winske 
[1987]  has  argued  that  the  instability  is  the  lower  hybrid  drift  instability  exhanced  in  the 
expanding  plasma  situation  by  the  deceleration  of  the  plasma.  Earlier,  Batchelor  and 
Davidson  [1975]  had  used  the  properties  of  the  instability  (claiming  g  ^  0)  to  explain  the 
flute  modes  observed  in  theta-pinch  compressions,  while  Okada  et  al.  [1981]  applied  similar 
theory  with  g  ^  0  to  interpret  their  laser  experiments.  Winske  [1987]  has  shown  that  the 
long  wavelength  limit  of  the  lower  hybrid  drift  instability  with  deceleration  [hereafter 


referred  to  as  the  LHDI]  reproduces  the  Hassam  and  Huba  [1987]  growth  rate  scaling.  The 
maximum  growth  rate  of  the  LHDI,  however,  occurs  at  very  short  wavelengths  ( kpe  ~  1, 
pe  =  electron  gyroradius),  which  does  not  seem  to  explain  the  longer  wavelengths  observed 
in  the  barium  release  or  the  laser  experiments.  Other  related  linear  theory  calculations 
and  simulations  [Galvez  et  al. ,  1987]  have  reached  the  same  conclusion,  although  their 
simulations  demonstrate  good  agreement  between  the  structure  wavelengths  observed  in 
their  simulations  and  the  wavelengths  at  maximum  growth  predicted  by  linear  theory. 
It  has  been  suggested  that  nonlocal,  nonlinear,  finite  beta,  or  velocity  shear  effects  may 
tend  to  inhibit  the  shortest,  fastest  growing  modes  from  dominating  at  late  times,  but  no 
conclusive  studies  have  yet  been  carried  out.  However,  until  one  understands  the  role  of 
these  short  wavelength  modes,  simulations  based  on  numerical  models  which  ignore  them 
e.g.,  hybrid  codes  with  massless  electrons,  or  single  fluid  codes]  are  suspect. 

In  our  previous  report  [Winske,  1987]  we  discussed  in  detail  the  linear  and  nonlinear 
theory  for  the  LHDI  and  presented  particle  simulations  of  dense  plasma  clouds  expanding 
into  a  tenuous  background  plasma.  While  a  number  of  cases  were  run  and  briefly  described, 
we  concentrated  there  on  one  particular  simulation  and  compared  in  detail  the  results  of 


the  calculations  with  theory.  Reasonably  good  agreement  was  obtained,  when  the  local 
parameters  were  well  known.  In  this  report  we  discuss  a  number  of  simulations  that  have 
been  done  and  emphasize  the  differences  in  the  results  as  various  parameters  are  changed. 
We  then  try  to  understand  the  behavior  of  the  instability  under  different  conditions  in 
terms  of  the  theory  of  the  LHDI.  In  conducting  such  studies  one  must  keep  in  mind  that 
the  parameters  are  not  all  independent.  For  example,  varying  the  ion  to  electron  mass  ratio 
affects  not  only  the  time  scale  for  the  instability  (~  ujJJj  ~  (me/m,) 1/2),  but  also  changes 
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the  magnetic  confinement  radius  ( Rb )  as  well  as  the  deceleration  rate.  In  addition  to 
changing  the  mass  ratio  (with  or  without  changing  we  have  also  varied  the  magnetic 
field  strength  (again  with  or  without  altering  Rb),  the  initial  velocity  profile,  and  the 
angular  spread  of  the  expanding  plasma.  Furthermore,  we  have  examined  the  direction  in 
which  the  flutes  bend  and  looked  for  their  bifurcation  at  late  times.  Overall,  the  results 
seem  to  be  in  good  qualitative,  if  not  quantitative,  agreement  with  the  theory.  We  have 
also  carried  out  a  rather  careful  analysis  of  the  simulations  to  see  whether  the  shortest 
waves  observed  in  the  simulations  are  physical  or  numerical  in  origin.  Unfortunately, 
we  find  that  electron  effects,  which  are  included  in  the  calculations  and  are  expected  to 
provide  a  minimum  scale  length,  do  not  prevent  even  shorter  wavelength  modes  from 
growing.  Explanations  based  on  both  physics  and  numerics  are  invoked  to  interpret  the 
results.  The  most  likely  explanation  is  that  another  instability,  the  electron  cyclotron  drift 
instability,  remains  unstable  under  most  conditions  of  the  simulation  and  is  responsible  for 
the  observed  short  wavelength  structures.  The  fact  that  the  observed  wavelengths  scale 
with  the  grid  size  explains  their  lack  of  magnetic  field  or  mass  dependence. 

The  outline  of  this  report  is  as  follows.  We  briefly  review  the  principal  results  of  the 
laser  experiments  of  Ripin  et  al.  [1987]  in  Section  2  and  of  the  theory  of  the  LHDI  in 
Section  3.  The  various  simulations  are  described  in  Section  4  and  interpreted  in  terms  of 
the  theory.  Section  5  addresses  numerical  issues,  and  Section  6  summarizes  the  results. 
The  overall  conclusion  is  that  the  shortest  wavelength  modes  which  are  determined  by  the 
numerical  grid  seem  unavoidable,  casting  doubts  on  any  scaling  studies  based  on  them. 
Given  such  restrictions,  the  best  course  of  action  seems  to  be  the  use  of  more  macroscopic 
codes  (e.g.,  hybrid  methods)  with  a  suitable  resistivity  algorithm,  or  some  other  numerical 
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trick,  to  try  to  flush  out  the  initial  perturbations  so  that  they  do  not  dominate  at  late 


times. 


2.  Review  of  the  Results  of  the  NRL  Laser  Experiment 

The  NRL  laser  experiments  were  conducted  with  a  Nd-laser  (l013Wy cm2, 30J)  onto 
a  small  A1  disk  target,  producing  a  plasma  mass  of  about  0.2  figm,  half  of  which  expanded 
forward  into  a  strong,  0. 1-1.0  T  magnetic  field.  The  outer  surface  of  the  expanding  plasma 
structures  into  field  aligned  (flute)  modes,  as  was  observed  by  means  of  a  gated  optical 
imager  camera.  Figure  1  (from  Ripin  et  ai.,  unpublished]  shows  an  example  of  the  phe¬ 
nomena.  This  particular  shot  was  taken  at  99ns  with  B=1T  at  somewhat  higher  laser 
energy  (367  J).  One  sees  that  the  plasma  has  expanded  from  the  target  (left  edge  of  the 
picture)  rather  symmetrically  and  the  outer  edge  of  the  plasma  has  structured,  produc¬ 
ing  finger-like  projections  which  are  fairly  regular  in  shape.  Other  measurements  indicate 
that  the  tips  of  the  flutes  are  freestreaming,  i.e.,  not  slowed  down  by  the  magnetic  field, 
although  at  late  times  the  flutes  acquire  a  curved  shape,  which  is  in  the  electron  E  x  B 
or  electron  cyclotron  sense  of  rotation.  These  is  also  a  suggestion  that  the  flutes  bifurcate 
in  time  to  maintain  roughly  constant  wavelength,  although  the  measurements  are  made 
at  different  times  on  different  shots.  The  wavelength  of  the  flute  modes  does  not  appear 
to  be  a  function  of  the  applied  magnetic  field.  This  is  in  contrast  to  the  earlier  Japanese 
experiments  [Okada  et  al.,  1981]  which  showed  that  the  wavenumber  of  the  unstable  modes 
scaled  with  the  magnetic  field,  ~  B0  8. 

It  should  be  noted  that  in  these  experiments  the  size  of  the  plasma  (~  Rb)  was 
comparable  to  the  ion  gyroradius,  and  much  larger  than  the  wavelength  of  the  observed 


structures.  For  example,  for  B=1T,  Rb  is  about  3.3cm  and  pi  was  estimated  to  be  about 
1.4cm,  while  at  B=0.1T,  Rb  —  15cm,  and  pt-  ~  14cm.  Similarly,  the  time  scale  for  the 
instability  (<  10  _8  sec)  was  much  less  than  the  ion  gyroperiod.  It  is  interesting  to  note  that 
the  barium  releases  in  the  magnetotail  (Bernhardt  et  al. ,  1987]  had  much  different  physical 
dimensions  and  values  for  the  parameters,  but  again  A  <  ft  ~  Rb-  In  both  situations, 
as  expected,  the  expanding  plasma  pushes  out  the  magnetic  field  to  form  a  diamagnetic 
cavity.  Magnetic  probe  measurements  on  the  laser  experiment  show  that  about  30%  of 
the  field  is  excluded  within  the  plasma  front.  In  both  experiments  high  frequency  electric 
field  noise  was  also  detected  (Ripin  et  al.,  1987;  Gurnett  et  al.,  1986].  Finally,  the  laser 
experiments  also  indicate  that  the  onset  of  the  instability  (in  time  and/or  distance  from 
the  target)  was  not  a  strong  function  of  the  applied  magnetic  field. 

3.  Review  of  Linear  Theory  of  the  LHDI 

A  detailed  discussion  of  the  lower  hybrid  drift  instability  in  its  generalized  form  is  given 
in  Winske  [1987],  where  the  linear  dispersion  equation  is  derived  and  numerous  analytic 
and  numerical  solutions  are  presented.  We  restrict  ourselves  here  to  a  few  relevant  features 
that  will  be  needed  to  interpret  the  simulations. 

The  physical  mechanism  which  gives  rise  to  the  instability  is  shown  in  Figure  2.  For 
times  short  compared  to  an  ion  cyclotron  period  (f ),t  <  1),  the  ions  can  be  taken  as 
unmagnetized.  As  they  expand  radially  outwards  across  the  magnetic  field,  the  ions  are 
constrained  by  the  electrons  which  are  tied  to  the  magnetic  field.  An  electric  field  is  then 
set  up  at  the  surface,  which  retards  the  ions  and  pulls  out  the  electrons.  In  the  radial 
direction,  the  expansion  thus  continues,  with  the  electrons  and  ions  being  compressed  into 


m 


a  thin  coupling  shell.  Azimuthally,  the  radial  electric  field  drives  the  electrons  via  E  x  B 


relative  to  the  ions,  giving  rise  to  the  instability. 


In  order  to  derive  a  linear  dispersion  equation  that  describes  the  local  properties  of  this 


azimuthally  directed  instability,  one  applies  the  response  of  the  unmagnetized  ions  and  the 


strongly  magnetized  electrons  to  perturbations  in  the  0  direction.  The  effect  of  the  radial 


expansion  on  the  local  dispersion  equation  enters  only  through  the  E  x  B  drift  velocity  of 


the  electrons.  Generally,  the  radial  electric  field  consists  of  two  parts.  One  piece  is  due 


to  the  slowing  of  the  plasma  as  it  pushes  against  the  magnetic  field  (characterized  by  a 


deceleration  g  =  —dVr/dt.  The  other  piece  is  due  to  the  density  gradient  (en  =  ~nldni/dr ), 


which  is  the  dominant  part  when  the  plasma  is  not  accelerating.  Assuming  for  simplicity 


cold  ( Te  =  0)  electrons,  the  radial  electric  field  can  be  written  as 


eEr  =  -rriig  -  Titn 


where  mt  =  ion  mass,  !T,=ion  temperature=|miV?,  so  that  the  E  x  B  velocity,  V#,  in  a 


magnetic  field  Ba  is 


Vp--c1l-  r!h£l  +  Iil.e  —  y  +  v 

E~  Bo  ~  eB0  +  eB0en~Va  +  Vn 


The  usual  lower  hybrid  drift  instability  corresponds  to  Vg  =  0;  the  plasma  expansion 


in  both  the  laser  and  the  simulations  to  be  described  later  is  characterized  by  Vg 


Vn.  In  electrostatic  simulations,  e.g.,  Galvez  et  al.  [1987],  the  part  of  the  electric  field 


characterized  by  g  is  not  a  consequence  of  compressing  the  magnetic  field,  but  rather 


results  from  charge  separation.  A  simple  model  for  g  in  this  situation  has  yet  to  be 


derived,  but  it  can  be  shown  from  the  simulations  that  Vg  ~  Vn. 


Proceeding  in  the  usual  fashion,  the  linear  dispersion  equation  can  be  derived  [Winske, 


S 


1987j.  It  is  not  written  down  here,  as  by  itself  it  is  not  very  revealing.  Analytical  approxi¬ 
mations  can  be  obtained  in  various  limits,  however,  which  are  more  useful  here.  Generally, 
unstable  solutions  occur  for  0  <  k  <  k±.  For  k  ~  0  (long  wavelength  limit),  the  growth 
rate  (-7)  scaling  can  be  shown  to  that  of  Hassam  and  Huba  [1987],  and  ur,  the  real  fre¬ 
quency,  is  very  small.  Maximum  growth  of  the  instability,  however,  occurs  at  much  shorter 
wavelength,  k  ~  km  <  k  1.  There  is  some  question  at  present  whether  the  modes  with  max¬ 
imum  growth  rate  are  relevant  to  either  the  laser  experiment  or  other  situations,  such  a s 
AMPTE.  They  are,  however,  excited  and  are  important  in  the  simulations,  if  sufficient  res¬ 
olution  exists  for  their  occurrence.  In  this  report  we  shall  be  primarily  concerned  with  the 
waves  of  maximum  growth  and  how  they  scale,  and  unless  otherwise  noted,  “by  comparing 
with  linear  theory”  we  will  mean  “by  comparing  with  linear  theory  for  the  fastest  growing 
waves”.  The  fastest  growing  modes  are  also  characterized  by  o/r  ~  7  ~  u jlh-  The  fact 
that  the  structures  observed  in  the  laser  experiments  (and  the  barium  releases)  show  little 
azimuthal  movement  suggests  u/r  ~  0.  This  might  indicate  that  the  flute  modes  seen  in 
the  experiments  have  become  nonlinear  with  properties  much  different  than  they  possesses 
at  earlier  time,  or  that  hardier,  low  frequency  (and  longer  wavelength)  modes  have  taken 
over  at  later  time.  In  the  simulations,  the  mode  number  is  roughly  constant  and  the  modes 
appear  to  have  a  larger  frequency  at  early  times,  suggesting  a  nonlinear  frequency  shift 
occurs  later,  as  happens  for  the  usual  lower  hybrid  drift  instability  without  deceleration 
VVinske  and  Liewer,  1987].  Nevertheless,  comparison  with  the  strongest  linearly  growing 
modes  seems  a  good  place  to  start. 


The  linear  dispersion  equation  yields  the  following  solutions  for  the  most  unstable 
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where  L,  =  c/wt,  H,  =  eB0/mtc,  and  =  Alfven  speed 
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where  u;r,//  =  wj/(  1  +  ^/n;!)1/2  is  the  lower  hybrid  frequency,  and  /?,•  =  8znTi/Bg.  For 
V3  =  0,  this  reduces  to  the  well  known  results  for  the  usual  lower  hybrid  drift  instability; 
when  Vg  7  0,  the  presence  of  the  deceleration  decreases  k  and  increases  wr  and  7. 

For  Fg  »  Vn,  Equations  (3)  and  (4)  reduce  to 
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agreement  with  (5)  and  (6)  and  at  least  qualitative  agreement  with  (7)  and  (8). 

While  these  results  are  interesting,  in  order  to  compare  later  with  the  simulations, 
they  need  to  be  expressed  somewhat  differently.  In  the  simulation  we  will  vary  :ue/fie 
and  mt/m,  and  see  how  the  mode  number  or  wavenumber  (i.e.,  k)  and  the  growth  time 
(~  7_1)  change.  Equations  (7)  and  (8)  then  give 
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where  ~  (Hen,)1/2  for  w2/fi2  »  1,  and  we  have  expressed  k  and  7  in  terms  of  we, 
the  plasma  frequency,  because  this  is  the  normalizing  quantity  in  the  simulations.  These 
results  demonstrate  a  clear  dependence  for  both  7  and  k  on  we/f2e  and  m,/me. 

However,  we  have  yet  to  relate  g  to  other  parameters.  This  is  more  important  for 
(9)  than  (10),  because  of  the  very  weak  dependence  on  Vg  in  (10).  The  deceleration  g 
can  be  related  to  other  parameters  assuming  radial  pressure  balance  [Huba,  1987].  In  two 
dimensions  this  can  be  expressed  as 
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where  N  is  the  number  of  cloud  ions  =  ncnr2L  ( r0  =  initial  cloud  radius,  nc=initial  cloud 
density,  L=  scale  length  in  third  (z)  direction).  Maximum  g  occurs  at  the  end  of  the 
expansion,  where  vr  =  0  and  r  =  Rq,  the  magnetic  confinement  radius: 
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where  jj2  =  47r nce2/mt.  Thus, 
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■  h ich  indicates  a  weaker  dependence  on  both  mxjmt  and  u/e/fie  than  is  implied  by  (9). 

The  difficulty  with  this  approach  is  that  the  instability  forms  at  early  times,  long 
efore  the  plasma  reaches  RB.  This  suggests  another  model,  which  leads  to  a  different 
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scaling.  Assume  at  early  times  that  the  ions  stream  ahead  of  the  electrons  a  small  distance 
Ar  and  exclude  the  magnetic  field.  Pressure  balance  is  maintained  by  the  electric  field  Er 
set  up  in  the  process  (ignoring  the  ion  pressure  contribution): 
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where  h  is  the  ion  charge  density  in  the  layer,  or 
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Since  (1)  implies  g  ~  — eEr/m.i ,  (9)  becomes 
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independent  of  both  we/ne  and  mi/me.  Thus,  various  scaling  for  the 
possible,  depending  on  the  assumed  physics  of  the  expansion. 


(16) 


(17) 


(18) 


(19) 

wavevector  are 


4.  Simulations 


3 


6 

x 

i"’ 

C 

\ 

| 

L‘ 


I 


This  section  consists  of  three  parts.  Part  A  is  an  overview  of  the  simulations,  reviewing 
previous  work  and  describing  the  methods.  In  Part  B  numerous  simulations  in  which  the 
initial  parameters  are  varied  are  discussed,  and  the  results  are  interpreted  in  terms  of 
the  theory  previously  presented.  Part  C  considers  other  issues  important  to  the  laser 
experiment,  such  as  the  late  time  behavior  of  the  flutes. 

A.  Overview 

The  simulations  are  similar  to  those  discussed  in  Winske  [1987!.  They  were  done  with 
the  fully  electromagnetic,  two-dimensional  particle  code  WAVE.  The  principal  difference  is 
that  here  there  is  no  (hot)  background  plasma.  A  dense  plasma  cloud  (typically  ncjn0  = 
25.  radius  r0  =  c  jje)  expands  at  VD  =  vA  {v\  =  3^/4™^,)  outward  against  a  uniform 
magnetic  field  B0,  which  is  in  the  z  direction  perpendicular  to  the  expansion.  Here  n0 
is  the  reference  density  that  defines  the  system  size,  typically  Lx  —  Ly  =  25 c/we,  where 
jjf  —  Airn0e2  ' mf,  with  the  magnetic  field  given  in  terms  of  f2e/we  =  eB0/mecu>e  (typically 
jje  =  0.2).  Usually  1252  cells  are  used  to  represent  the  system,  and  40000  ions  and 
electrons  describe  the  expanding  plasma.  Nominal  parameters  include:  m,/me  =  100, 
ty  c  =  0.001,  and  vtjc  =  0.01,  with  VD  =  0.02c.  The  choice  of  these  values  represents  a 
compromise:  one  must  keep  the  system  size  reasonable,  yet  resolve  scales  of  interest.  We 
thus  tolerate  some  numerical  heating  of  the  electrons  early  in  the  simulation  (since  < 
A):  previously  this  was  minimized  with  the  addition  of  a  hot  background  component. 
Elimination  of  the  background  does  not  change  the  results  significantly,  as  will  be  shown 
later,  and  speeds  up  the  calculations  considerably.  Finally,  we  note  that  in  contrast  to  the 
calculations  here,  the  simulations  of  Galvez  e£  al.  [1987]  are  electrostatic  and  involve  the 
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thermal  expansion  of  hot,  large  (many  ion  gyroradii)  plasmas. 


An  example  of  one  such  simulation,  with  the  nominal  parameters  given  above,  is 
shown  in  Figure  3.  Displayed  are  the  ions  in  real  space  (x-y)  at  various  times  [in  the  figures 
-^LHt  =  (nen,)1/2t].  At  early  time  the  outer  edge  of  the  expanding  cloud  is  slowed  by  the 
magnetic  field,  while  the  inner  edge  is  not,  compressing  the  cloud  into  a  thin  “coupling 
shell”.  At  somewhat  later  time,  modes  appear  on  the  surface.  The  waves  continue  to  grow 
as  the  plasma  expands  and  eventually  the  plasma  diffuses  through  the  magnetic  field. 
Compared  to  the  case  with  a  background  plasma  present  [Figure  11,  Winske  (1987)],  there 
is  little  difference  in  the  size  or  shape  of  the  expanding  cloud.  The  modes  in  the  present 
simulation  are  somewhat  more  diffuse  and  of  slightly  shorter  wavelength  (higher  mode 
number  m),  but  overall  the  differences  are  minor. 

B.  Parameter  Variation 

We  next  compare  the  results  of  various  simulations  in  which  some  of  the  input  pa¬ 
rameters  are  varied.  Without  doing  detailed  comparisons  of  each  case  with  linear  theory 
(which  requires  accurate  local  measurements  of  quantities  such  as  the  relative  electron-ion 
drift  speed),  we  will  attempt  to  interpret  the  results  in  terms  of  theory  and  relate,  where 
possible,  to  the  experimental  results.  Generally,  the  parameters  will  be  varied  about  the 
“Standard  Case”  described  in  the  previous  subsection  and  snapshots  of  the  ion  density 
contours  at  a  particular  time  rather  than  a  sequence  of  plots  will  be  shown. 

The  first  case  considers  variation  with  ion  mass,  m,/me,  the  other  parameters  kept 
constant.  The  results  are  shown  in  Figure  4:  top  panel,  m,/me  =  25  at  wet  =  175;  middle 
panel,  standard  case,  m,/me  =  100,  uet  —  225;  bottom  panel,  m,/me  =  400,  at  uet  =  300. 
The  times  are  chosen  so  that  the  instability  is  at  about  the  same  stage  of  development  in 


each  case.  It  is  clear  that  the  heavier  ion  mass  implies  that  the  instability  develops  more 
slowly  (in  terms  of  clock  time,  i.e.,  u>et)  so  that  the  size  of  the  plasma  shell  is  larger  and 
the  mode  number  m  is  larger.  However,  the  wavenumber  in  the  three  cases  is  not  much 
different,  ck/uie{=  2n/X)  =  4.9, 3.9, 4.0,  respectively  for  m,/me  =  25,100,400,  and  given 
the  difficulty  in  counting  the  modes,  the  results  suggest  that  there  is  no  mass  dependence 
on  the  wavenumber.  Note  that  expressing  the  time  in  terms  of  uLH  [~  (Hen,)1/2],  the 
three  frames  in  Figure  4  correspond  to  unit  =  7, 4.5, 3.  [One  cannot  say,  however,  how  the 
actual  growth  rates  vary  in  the  simulations  by  just  looking  at  these  plots,  since  the  time 
over  which  there  is  strong  exponential  growth  of  the  waves  is  very  limited,  e.g.,  see  Figure 
15  in  Winske  (1987)].  As  noted  in  the  previous  section,  the  wavenumber  may  have  a  rather 
strong  [mi/ me)1/2  dependence  on  mass,  moderate  dependence  (m./me)1/4,  or  no  depen¬ 
dence,  depending  on  conditions  and  the  physics  model  assumed.  Thus,  the  simulations 
are  consistent  with  theory,  assuming  (16)-(19)  is  correct.  An  alternative  possibility  to  be 
kept  in  mind,  which  is  discussed  in  detail  in  the  next  section,  is  that  the  mode  structure 
is  instead  limited  by  the  spatial  resolution  of  the  grid;  i.e.,  the  maximum  wavelength  is 
k  ~  1/Ax  ~  5 we/c. 

In  similar  fashion,  results  showing  the  variation  of  Cl e/we  are  displayed  in  Figure 
5:  top  panel,  He/u;e  =  0.1  at  wef  =  300;  middle  panel,  O e/we  =  0.2  at  =  225, 
the  standard  case;  bottom  panel,  Cl e/we  =  0.4  at  ojet  =  175.  Again,  the  ion  density 
contours  are  shown  at  a  comparable  stage  of  development  of  the  instability.  As  before, 
in  terms  of  normalized  time,  the  times  are  reversed,  with  the  top  panel  actually  being 
earliest  {uLHt  —  3)  compared  to  the  other  cases  (u )LHt  =  4.5  and  7,  respectively).  At 
larger  values  of  Cle/oje  the  instability  develops  faster,  hence  the  plasma  has  a  smaller 


radius  and  mode  number.  Again,  kc/uit  =  3.8, 4.0, 4.4  for  the  cases  ne/we  =  0.1, 0.2, 0.4, 
respectively.  As  with  the  mass  dependence,  linear  theory  for  Cle/u>e  variations  predicts 
difference  dependences  in  different  regimes.  Again,  the  estimate  of  k  based  on  assuming  a 
charge  separation  layer  (16-19)  gives  no  fie/we  dependence.  In  contrast,  in  the  simulations 
of  Galvez  et  al.  [1987],  which  have  g  ~  0,  the  wavenumbers  scale  like  Eq.  (5).  One 
recalls  from  the  experimental  results,  that  there  is  also  some  ambiguity.  Ripin  et  al.  [1987] 
indicate  no  variation  of  wavelength  with  B  (i.e.,  Uejut),  consistent  with  Eqs.  (16)- (19), 
while  Okada  et  al.,  [1981]  report  k  ~  B0-8,  consistent  with  (11)-(15).  Again  the  differences 
may  suggest  the  two  experiments  are  run  in  different  regimes.  For  example,  the  modes 
in  Okada’s  experiment  occur  at  maximum  expansion  [consistent  with  ( 1 1-15)],  while  in 
Ripin’s  experiment,  the  structures  appear  earlier,  more  like  (16)-(19). 

In  order  to  further  check  the  scaling  and  to  normalize  the  differences,  we  have  also 
done  simulations  in  which  Rb  is  kept  constant.  For  example,  in  Figure  6  we  have  varied 
both  m,/me  and  O e/u;e  together:  top  panel,  ne/we  =  0.1,  m,/me  =  25;  middle  panel, 
ne/we  =  0.2,  m,/me  =  100,  standard  case;  bottom  panel,  fie/u/e  =  0.4,  mi/me  =  400. 
In  each  case  um  is  identical,  with  the  results  at  uet  =  300,  uujt  =  6.  In  this  case,  the 
size  of  the  clouds  are  comparable  (~  i?s)>  the  instability  is  about  in  the  same  state  of 
development,  and  the  mode  numbers  are  roughly  the  same.  This  result  is  consistent  with 
both  the  scaling  (15)  and  (19). 

Another  way  to  compare  the  simulations  is  to  vary  fle/we  and  Vp/c  to  again  keep  Rb 
constant.  Figure  7  displays  the  ion  density  contours  for  three  cases  all  at  =  6:  top 

panel,  fie/we  =  0.1,  VD/c  =  0.01;  middle  panel,  ne/we  =  0.2,  V^/c  =  0.02,  standard  case; 
bottom  panel,  Qe/uie  =  0.4 ,Vd/c  =  0.04.  Again,  the  instability  has  developed  to  about 


the  same  point  (the  top  panel  looks  weaker,  but  may  indicate  a  lower  saturation  level),  the 
size  of  the  plasmas,  and  the  wavelengths  are  comparable  in  all  cases.  Again,  the  results 
are  consistent  with  both  (15)  and  (19). 

The  effect  of  varying  the  initial  conditions  has  also  been  investigated.  Figure  8  shows 
the  standard  case  in  which  all  the  particles  have  the  same  initial  radial  velocity  VD  and  a 
second  case  where  a  self-similar  initial  radial  velocity  (\/2 VDr/r0)  distribution  was  used. 
The  y/2  is  included  so  that  the  total  kinetic  energy  is  the  same  in  each  distribution. j  The 
figure  shows  that  both  cases  give  similar  results,  about  the  same  mode  number  and  level 
of  charge  bunching.  Because  of  the  velocity  spread,  however,  some  ions  remain  behind  in 
the  expansion  in  the  self-similar  case,  and  at  later  times,  it  spreads  to  larger  radius. 

The  initial  angular  spread  of  the  expansion  has  also  been  varied.  Figure  9  displays 
three  self-similar  initial  distributions  ( Vr  —  VDr/r0)  with  angular  spread  ±0o  (180°,  90°,  45°) 
and  one  VT  —  VD  case  with  60  —  90°,  all  at  wet  —  300.  The  three  self-similar  cases  show 
some  variation  with  60:  large  0O  gives  a  slightly  longer  wavelength  and  the  flutes  are  some¬ 
what  more  elongated.  Both  effects  could  be  caused  by  the  larger  deceleration  expected 
when  60  is  larger.  Also  note  in  the  partial  shells  that  the  instability  is  also  present  on  the 
rear  side  of  the  expansion  at  about  the  same  wavelength.  In  this  case  the  density  gradient 
drift  is  very  small,  as  the  ion  thermal  speed  is  negligible,  and  the  deceleration  drift  is  about 
the  same  as  on  both  sides  of  the  shell.  Finally,  the  bottom  panel  shows  the  constant  initial 
velocity  case  with  0o  =  90°.  Three  effects  are  interesting:  (1)  there  is  little  slowing  of  the 
shell,  (2)  the  waves  that  form  on  the  surface  are  quite  small,  and  (3)  the  entire  shell  has 
rotated  slightly  in  the  counterclockwise  direction.  Because  the  ions  are  unmagnetized,  this 
rotation  must  result  from  the  ions  following  the  electrons,  which  gives  the  correct  sense. 


The  fact  that  there  is  little  deceleration  probably  results  from  the  fact  that  little  B  field 
get  excluded;  the  flux  can  just  flow  around  the  shell.  The  small  deceleration  could  imply 
smaller  growth  of  the  instability  and  lower  saturation  level  [cf.,  Winske  (1987)  for  the  case 
of  expansion  in  a  slab  geometry]. 

C.  Other  Issues 

Besides  demonstrating  the  growth  of  the  instability,  the  laser  experiments  of  Ripin 
et  al.  1987]  exhibited  very  interesting  late  time  behavior,  including  the  bifurcation  of  the 
flutes  and  their  bending  in  time.  The  issue  of  whether  the  freestreaming  flute  tips  can 
break  off  and  continue  propagating  across  the  magnetic  field  is  also  important  for  the 
VHANE  situation.  Generally,  the  simulations  [Winske,  1987]  show  that  by  the  time  the 
flutes  are  visible  in  real  space,  wave  growth  has  saturated.  Thus,  much  of  what  is  observed 
is  nonlinear  phenomena  that  is  not  described  by  the  theory  presented  earlier.  More  work, 
both  experimental  and  theoretical,  will  be  needed  to  understand  the  observations.  The 
results  discussed  here  are  included  to  demonstrate  that  the  particle  simulations  can  be  use¬ 
ful  for  these  studies,  although  hybrid  (massless,  fluid  electrons)  [e.g.,  Brecht  and  Thomas, 
unpublished]  will  be  better  because  of  the  larger  spatial  and  temporal  scales  that  can  be 
used. 

The  top  panel  of  Figure  10  shows  the  ion  density  contours  for  a  run  presented  earlier 
[Figure  5,  bottom  panel]  at  later  time.  The  Y-shaped  structure  at  the  top  of  the  figure 
suggests  a  bifurcated  flute  mode.  The  lower  panels  show  snapshots  of  the  density  at  earlier 
times,  in  reverse  order.  It  is  clear  from  these  pictures  that  the  flute  tip  does  not  split,  but 
that  the  modes  just  evolve  in  time,  and  the  mode  number  remains  pretty  much  unchanged. 
Similarly,  comparing  the  middle  panel  of  Figure  4  with  the  top  panel  of  Figure  8  (i.e.,  the 
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standard  case  at  two  times)  basically  shows  that  the  mode  number  does  not  change  very 
much.  Sometimes,  e.g.  the  bottom  panel  of  Figure  7,  the  contour  plotter  groups  together 
two  smaller,  separate  peaks  into  one  larger  structure.  Overall,  the  constancy  of  the  mode 
number  is  consistent  with  older  simulations  of  the  lower  hybrid  drift  instability  in  slab 
geometry  i  VVinske  and  Liewer,  1978]  as  well  as  recent  mode  coupling  theory  [Drake  et  al., 
1984  .  A  difference  here,  however,  is  the  continuous  expansion  of  the  plasma  allows  the 
waves  to  be  driven  to  higher  amplitudes  where  stochastic-like  decoupling  of  the  electrons 
from  the  magnetic  field  [Drake  and  Lee,  1981]  occurs.  Thus,  the  electron  behavior  of  the 
usual  mode  coupling  model  that  involves  V B  resonaces  may  be  expected  to  be  modified. 


Another  interesting  feature  of  Figure  10  is  the  direction  that  the  flutes  bend.  Except 
for  the  Y-shaped  structures,  the  arms  generally  spiral  clockwise,  i.e.,  in  the  ion  sense  of  ro¬ 
tation  (recall  the  B  field  is  outward).  This  is  characteristic  of  most  of  the  simulations,  e.g., 
Figure  3.  The  behavior  is  not  due  to  ion  gyromotion,  because  the  ions  are  unmagnetized 
in  the  simulations.  From  examining  the  temporal  sequence  in  Figure  10,  it  is  evident  that 
these  curved  arms  are  all  rotating  at  a  slow  rate  in  the  electron  (counterclockwise)  sense. 
One  expects  the  arms  to  bend  clockwise  if  the  ends  of  the  arms  are  just  moving  outward, 
but  the  bases  rotate  counterclockwise  as  the  plasma  there  is  carried  by  the  waves  which 
have  phase  velocities  in  the  E  x  B  (counterclockwise)  direction.  At  early  times,  ur  ~  ujlh 
but  the  waves  do  not  affect  the  particles  yet,  while  at  late  times  the  waves  experience  a 
nonlinear  phase  shift  so  that  wr  *C  ^lh- 

Although  this  is  a  reasonable  interpretation,  Figure  11  shows  another  case  where  the 
arms  spiral  in  the  electron  sense  of  rotation.  This  case  is  characterized  by  a  very  rapid 


growth  of  the  instability  (He/u;e  =  0.4  and  nc/n0  =  100),  which  clumps  the  density  at  very 
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early  times  and  then  a  continued,  unconfined  expansion  across  the  magnetic  field  ensues. 

Finally,  Figure  12  shows  another  dense  plasma  ( ncjn0  =  100,  r0  =  2 c/we)  expansion 
into  a  weak  (fie/u;e  =  0.1)  field.  In  this  case  strong  density  bunches  are  formed  at  early 
time.  Again,  here  there  is  not  much  confinement  of  the  plasma.  The  rapid  growth  of 
the  plasma  waves  to  large  amplitude  allows  the  plasma  to  stream  through  the  field.  The 
interesting  aspect  of  this  example  is  that  the  density  clumps,  once  formed,  do  not  diffuse. 
Instead,  the  plasma  propagates  across  the  field  as  discrete  plasmoids.  While  one  has  to 
be  cautious  in  interpreting  the  simulation  results,  because  artificial  parameters  (mt/me, 
ne/a;e,  etc)  allow  the  wave  turbulence  and  hence  the  anomalous  resistivity  (diffusion)  to  be 
artifically  large,  it  does  suggest  that  interesting  behavior  can  occur  under  some  conditions. 

5.  Analysis 

In  the  previous  section  it  was  suggested  that  the  lack  of  ion  mass  or  magnetic  field 
strength  dependence  of  the  wavenumber  of  the  observed  flute  modes  was  not  physical 
and  instead  related  to  numerical  issues.  Here  we  examine  this  question  in  some  detail 
and  conclude  that  the  modes  in  the  simulation  are  always  determined  by  the  grid  size. 
However,  the  modes  are  physical  and  are  likely  due  to  a  shorter  wavelength  instability,  the 
electron  cyclotron  drift  instability.  The  consequences  of  this  result  are  discussed  in  the 
conclusions. 

Previous  simulations  of  plasma  expanding  across  magnetic  fields  in  which  electron 
effects  were  ignored  [Shonk  and  Morse,  1968;  Brecht  and  Thomas,  unpublished,  1987; 
Huba  et  al. ,  1987]  all  indicated  that  the  shortest  allowable  modes  grow  fastest.  This  is 
consistent  with  linear  theory  for  the  lower  hybrid  drift  instability  (  or  equivalently  the 


theory  for  the  unmagnetized  ion  Rayleigh  Taylor  instability  [Hassam  and  Huba,  1987])  in 
that  at  long  wavelengths,  -7  ~  k.  However,  when  electrons  are  included  in  the  analysis,  a 
sharp  maximum  in  the  growth  rate  curve,  given  by  (3),  occurs;  at  larger  k,  growth  falls 
rapidly  to  zero.  One  thus  expects  that  simulations  which  include  full  electron  dynamics 
would  show  that  modes  at  wavelengths  shorter  than  this  cutoff  to  be  suppressed.  This 
can  be  easily  checked  by  varying  the  grid  spacing  in  the  simulation.  Figure  13  shows  the 
results  of  this  test  for  the  Standard  Case  of  Figure  3  for  three  different  grid  spacings: 
(A)  A  =  0.125 c/we,  (B)  A  =  0.0625,  and  (C)  A  =  0.03125  (in  Figure  3,  A  =  0.2).  In 
each  case  we  have  kept  the  number  of  simulation  particles  the  same  and  varied  only  the 
size  of  the  simulation  region;  thus  cases  with  smaller  A  were  run  a  shorter  time.  The 
figure  demonstrates  clearly  that  more  modes  appear  a s  the  cell  size  shrinks  and  that  the 
shortest  modes  continue  to  dominate.  Also,  the  waves  tend  to  appear  at  earlier  times. 
Unfortunately,  we  are  unable  to  continue  the  fine  mesh  runs  in  time  without  rezoning, 
which  would  then  suppress  the  shortest  wavelengths. 

Several  explanations  for  this  behavior  are  possible.  The  first  possibility  is  that  the 
plasma  properties  change  with  grid  size.  For  example,  the  magnitude  of  the  radial  electric 
field,  and  hence  VE,  could  depend  on  A  (i.e.,  A  ~  Ar  in  (19)];  the  unstable  waves  then 
correspond  to  the  same  instability,  but  under  different  conditions.  In  order  to  check  this, 
we  have  examined  the  magnitude  of  the  cross  field  drift  and  the  properties  of  the  modes 
expected  from  linear  theory.  The  results  are  presented  in  Table  1.  Generally,  the  results 
show  that  the  cross-field  drift  speed  is  about  the  same  in  each  case  at  the  same  time. 
Here  VE  is  determined  from  both  narrow  slices  in  x  and  y  of  the  electron  jx  and  jy  and 


the  density,  as  well  as  from  Vg  versus  6  phase  space  for  the  electrons.  This  result  rules 
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out  an  electric  field  that  depends  on  A  .  The  magnitude  of  the  electric  field  fluctuations 
increases  as  A  decreases,  but  that  it  likely  due  to  the  fact  that  the  number  of  particles 
per  cell  is  decreasing.  The  radial  electric  field  also  increases  with  time.  The  magnitude 
of  Er  and  its  variation  with  time  is  consistent  with  the  simple  model  of  the  exclusion  of 
the  magnetic  field  iHuba,  1987],  as  shown  in  our  earlier  study  [Winske,  1987].  Thus,  for 
each  run  the  magnitude  of  Vg  increases  with  time,  which  according  to  linear  theory,  means 
that  k  is  decreased.  However,  the  plasma  radius  increases  with  time  so  that  the  mode 
number  m=kr  stays  roughly  constant,  consistent  with  the  simulations.  For  the  case  with 
the  largest  cell  size,  the  dominant  mode  predicted  by  linear  theory  is  roughly  the  same 
as  that  assuming  the  shortest  modes  dominate  (i.e.,  km  ~  l/A),  but  for  the  smaller  A 
cases  there  is  a  clear  separation  between  the  mode  number  predicted  by  theory  and  that 
observed  in  the  simulations. 

The  discrepancy  between  simulation  and  theory  can  be  due  to  either  numerical  or 
physical  reasons.  So  a  second  explanation  is  that  the  results  are  merely  numerical,  due 
to  trying  to  simulate  an  initially  round  plasma  on  a  square  grid.  Because  the  two  can 
never  match  completely,  there  will  always  be  short  wavelength  ( k  ~  1/A)  perturbations 
Brecht,  private  communication  .  A  counter  argument  to  this  is  that  such  perturbations 
are  expected  to  be  strongly  damped,  according  to  the  linear  theory.  The  resolution  of 
this  issue  awaits  the  development  of  an  r -9  PIC  code  with  variable  spacing  in  the  radial 
direction,  which  we  hope  to  have  running  sometime  next  year. 

A  third  explanation  based  on  physical  considerations  primarily  concerns  whether  there 
are  shorter  wavelength  modes  which  can  grow.  A  likely  candidate  is  the  electron  cyclotron 
drift  instability  ECDIl  Forslund  et  al.,  19701,  which  is  another  mode  driven  by  the  cross 


field  drift  of  the  electrons  relative  to  the  ions.  This  instability  occurs  at  higher  frequency 
that  the  lower  hybrid  drift  instability  and  involves  the  interaction  of  ion  waves  with  electron 
Bernstein  modes,  i.e.,  u>  ~  nfle,n  =  1,2,3....  The  instability  is  known  to  be  strong  when 
VE  ~  ve,  the  situation  in  the  simulation;  moreover,  it  is  characterized  by  short,  wavelengths 
kpe  3>  1.  In  Figure  14  we  display  solutions  of  the  linear  dispersion  relation  for  the  LHDI 
(Tj  =  0)  and  the  ECDI  for  the  simulation  parameters  of  Table  1  at  o oet  =  25.  Only  the 
two  lowest  branches  of  the  ECDI  are  shown.  Maximum  growth  (7m  ~  0.047we)  occurs  at 
much  larger  k  ( kc/ue  ~  160),  with  the  instability  persisting  for  even  much  larger  values 
of  k  (higher  n).  At  higher  mode  numbers  the  qr  verus  k  curves  broaden  so  that  eventually 
the  unstable  bands  overlap.  The  large  gap  between  the  LHDI  and  the  n=l  branch  of 
the  ECDI  decreases  when  finite  electron  temperature  effects  are  included  in  the  LHDI 
dispersion  relation.  The  conclusion  fron  this  graph  is  that  there  are  real,  short  wavelength 
unstable  linear  modes  which  can  grow. 

The  principal  difficulty  is  whether  such  unstable  waves  can  grow  to  significant  ampli¬ 
tude.  Generally,  the  ECDI  is  dismissed  as  being  unimportant  [Lampe  et  al.,  1972]  because 
it  saturates  at  low  amplitudes.  This  is  because  the  instability  depends  on  the  cyclotron 
resonances  of  the  electrons  and  magnetic  field  fluctuations  or  gradients  tend  to  destroy  the 
resonances  whenever  the  instability  grows  to  even  small  amplitudes.  However,  because  of 
the  persistence  of  the  electric  field  in  time,  the  instability  may  be  driven  to  larger  ampli¬ 
tudes  than  expected,  just  as  the  LHDI  is  driven  to  levels  far  exceeding  the  ion  trapping 
limit  usually  expected,  in  this  expanding  geometry  [Winske,  1987],  The  ECDI  (and  the 
LHDI  to  a  lesser  extent)  is  suppressed  when  0  [Lashmore-Davies  and  Martin,  1973]. 
Simulations  with  the  magnetic  field  tilted  so  that  Bx  0,  however,  show  no  real  differ- 


ence.  This  lack  of  stabilization  may  be  an  artifact  of  the  two-dimensional  nature  of  the 


simulation,  similar  to  the  observed  lack  of  shear  stabilization  of  the  LHDI.  Alternatively, 
the  finite  k ,  stabilization  may  be  removed  when  electromagnetic  effects  are  included,  as 
happens  when  the  instability  involves  ion  resonances  [Brecht,  1982], 

Finally,  one  asks  whether  there  are  regimes  where  the  shortest  wavelengths  can  be 
resolved.  Clearly,  if  one  starts  with  large  plasmas  (r0  ;»  p,  ~  Rb),  the  most  unstable 
modes  are  well  resolved  [Galvez  et  al.,  1987],  However,  for  the  laser  case,  one  wants 
to  study  the  opposite  situation,  where  r0  p  ~  Rb-  While  one  can  have  reasonable 
parameters  in  terms  of  ion  quantities,  e.g.,  Vp  ~  r0  -C  pt,  because  of  the  compression 
of  scales  by  choosing  uie/Cle  and  m,/me  small,  one  is  always  in  a  regime  where  Vg  is  not 
negligible  compared  to  ve  and  the  electron  cyclotron  drift  instability  always  seems  to  occur. 
For  example,  increasing  Te/T,  while  reducing  Vjs;/ve  does  not  improve  the  situation  very- 
much. 

6.  Conclusions 

In  summary,  in  this  report  we  have  described  results  from  a  number  of  particle  sim¬ 
ulations  of  dense  plasmas  expanding  at  the  Alfven  speed  into  a  vacuum  across  a  strong 
magnetic  field.  This  type  of  configuration  is  subject  to  a  strong  instability  that  produces 
flute  modes  on  the  outer  surface  of  the  expanding  plasma.  Such  effects  have  been  observed 
and  investigated  in  laser  produced  plasmas  [Ripin  ef  al.,  1987]  and  in  barium  releases  in 
space  [Bernhardt  ef  al.,  1987],  and  have  been  the  subject  of  much  recent  theoretical  work 
Hassam  and  Huba,  1987;  Huba  et  al.,  1987;  Winske,  1987;  Galvez  et  al.,  1987],  Our 
previous  report  [Winske,  1987]  suggested  that  the  lower  hybrid  drift  instability  enhanced 
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by  the  deceleration  of  the  plasma  across  the  magnetic  field  was  a  good  candidate  to  ex¬ 
plain  the  observations,  although  there  was  considerable  disagreement  between  the  observed 
wavelengths  and  those  corresponding  to  the  fastest  growing  waves  of  linear  theory. 

Our  purpose  here  has  been  to  describe  numerous  simulations,  varying  the  initial  con¬ 
ditions  in  a  systematic  fashion  to  understand  the  scaling  of  the  instability  with  various 
parameters  and  to  compare  with  theoretical  predictions.  Because  of  the  interconnection  of 
the  input  conditions  through  the  deceleration  process,  such  scaling  studies  are  not  trivial. 
We  have  varied  both  the  mass  ratio  and  the  magnetic  field  strength  (Qe/we)  indepen¬ 
dently.  in  which  case  the  size  of  the  plasma  varies,  as  well  as  with  other  parameters  so 
that  Rb  ^constant.  It  has  been  shown  that  the  wavelength  of  the  instability  is  indepen¬ 
dent  of  these  two  quantities,  contrary  to  the  standard  theory  in  which  the  deceleration  is 
related  to  the  compression  of  the  magnetic  field.  A  closer  examination  of  the  results  and 
the  simulation  methods  suggest  that  the  magnitude  and  time  variation  of  the  electric  field 
based  on  the  magnetic  field  compression  model  is  basically  correct.  The  wavelength  of  the 
instability,  however,  is  always  determined  by  the  grid  size.  For  wavelengths  shorter  than 
those  predicted  by  the  lower  hybrid  drift  instability,  unstable  modes  that  can  be  identified 
with  the  electron  cyclotron  drift  instability  persist. 

Given  the  limitation  of  the  minimum  wavelength  due  to  grid  effects,  we  have  neverthe¬ 
less  also  shown  some  variation  in  the  character  of  the  instability  with  the  initial  distribution 
in  radial  velocity  and  angular  spread.  We  have  seen,  for  example,  that  partial  shells  with 
self-similar  velocity  distributions  structure  on  the  inner  edge  with  about  the  same  prop 
erties  as  on  the  outer  edge,  but  partial  shells  with  constant  intial  velocity  distributions 
do  not  seem  to  structure  at  all.  In  addition,  a  number  of  late  time  (probably  nonlinear) 


effects  have  been  observed  in  the  simulations  that  are  relevant  to  the  experiment.  First,  we 
have  shown  that  the  tips  of  the  flute  modes  which  appear  to  bifurcate  are  actually  closely 
spaced  waves  that  are  diverging.  Overall,  there  seems  to  be  little  tendency  for  the  waves 
to  cascade  to  longer  or  shorter  wavelengths  in  the  simulations,  at  least  on  the  time  scale 
studied.  Exceptions  to  this  rule  are  cases  where  very  strong  instabilities  are  excited  and 
dynamic,  nonlinear  processes  (bending  or  breaking  off  of  the  flute  tips)  are  seen.  We  have 
observed  flute  tips  that  bend  in  the  electron  sense  of  rotation  in  the  magnetic  field  (like 
the  experiment)  as  well  as  in  the  ion  sense. 

The  work  included  in  this  report  constitutes  a  preliminary  study  of  the  basic  physics 
that  should  occur  in  the  laser  experiment.  So  far,  it  has  been  primarily  a  test  of  the 
capabilities  and  the  limitations  of  the  simulations  and  the  linear  theory  to  explain  them. 
And  it  has  pointed  out  areas  of  nonlinear  physics  that  need  further  study  to  explain 
the  data.  Furthermore,  it  has  demonstrated  the  difficulties  of  improving  the  numerical 
modeling  of  both  the  laser  configuration  and  the  AMPTE  barium  releases. 

The  overall  conclusion  of  this  study  points  to  a  fundamental  limitation  in  using  particle 
methods  to  model  the  laser  experiments.  The  original  purpose  of  employing  a  simulation 
model  with  full  electron  and  ion  dynamics,  which  sacrifices  some  reality  by  compressing 
ion  and  electron  scales,  was  to  estimate  the  importance  of  the  electrons  is  limiting  short 
wavelength  fluctations.  However,  a  detailed  analysis  of  the  simulations  demonstrate  that 
beyond  the  shortest  wavelengths  allowed  by  the  lower  hybrid  drift  instability,  other  unsta¬ 
ble  modes  due  to  the  electron  cyclotron  drift  instability  persist  down  to  the  shortest  scales 
realizeable.  The  simulations  thus  suggest  that  the  shortest  wavelengths  in  any  reasonable 
sized  simulation  will  be  dominated  by  numerical  grid  effects.  On  the  positive  side,  the 
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simulations  do  show  that  there  are  short  wavelength  modes  present,  even  at  very  early 


times,  which  heat  the  plasma  and  diffuse  it  through  the  magnetic  field.  The  radial  electric 
field  that  is  set  up  and  tends  to  increase  in  time  has  been  shown  to  be  consistent  with 
the  simple  model  of  magnetic  field  exclusion,  and  is  itself  not  dependent  on  the  grid  size. 
Because  of  the  continued  deceleration  of  the  plasma  by  the  magnetic  field,  the  azimuthai 
current  is  not  able  to  relax,  and  hence  the  instability  can  be  driven  to  rather  high  levels. 
The  net  effect  is  that  one  gets  a  rather  broad  (c/u>t),  turbulent  boundary  layer  where  the 
initial  waves  which  brought  this  about  should  be  washed  away.  Longer  wavelength,  slower 
growing  modes  may  then  develop  at  the  surface  of  this  broadened  layer. 

Given  the  fact  that  the  dominant  modes  in  the  simulation  at  early  times  at  determined 
numerically  (which  by  judicious  choice  of  the  grid  can  be  chosen  to  be  those  expected  from 
linear  theory),  one  would  like  to  run  the  simulation  long  enough  until  they  are  no  longer 
dominant.  With  the  full  particle  code  that  follows  electron  scales  that  is  generally  not 
possible,  although  some  of  the  runs  with  a  very  violent  instability  initially  (e.g.,  Figures  11 
and  12)  show  this  behavior.  A  better  approach  would  be  hybrid  calculations  [e.g.,  Brecht 
and  Thomas,  unpublished]  where  by  means  of  a  time  varying  resistivity  it  might  be  possible 
to  flush  out  the  initial  perturbations  and  then  let  longer  wavelength  modes  develop.  This 
washing  out  of  the  system  at  early  times  apparently  happens  in  the  experiments  and  may 
explain  why  the  observed  structures  seem  independent  of  the  magnetic  field  and  exhibit 
peculiar  nonlinear  behavior. 
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Table  1.  Simulation  parameters  for  runs  with  various  cell  sizes. 


(A) 

(B) 

(C) 

A(c/we) 

0.125 

0.0625 

0.03125 

uet 

25 

50 

25 

50 

25 

50 

nev$/n0c 

0.22 

0.25 

0.23 

0.25 

0.30 

0.25 

ne/n0 

15 

10 

15 

10 

16 

10 

v$/c 

0.015 

0.025 

0.015 

0.025 

0.017 

0.025 

^(phase  sp.) 

0.015 

0.025 

0.2 

0.025 

0.2 

0.025 

Vg/vA 

1.6 

2.4 

1.6 

2.4 

1.7 

2.4 

kLHDc/ue 

10.5 

8.0 

10.5 

8.0 

10.5 

8.0 

r(c/we) 

1.5 

2.0 

1.5 

2.0 

1.5 

2.0 

mLHD  =  kr 

16 

16 

16 

16 

16 

16 

^num  r  /  A 

12 

16 

24 

32 

48 

64 

Other  parameters:  Qe/we  =  0.2,  Te/T,  =  1,  m,/me  =  100,  ve/c  =  0.01,  and  c„c/u> 


Figure  Captions 

Figure  1.  Example  of  plasma  structuring  observed  in  the  NRL  laser  [ Ripin  et  al. ,  unpub¬ 
lished]. 

Figure  2.  Instability  mechanism. 

Figure  3.  Simulation  of  the  instability,  showing  the  ions  in  real  space  at  various  times. 

Figure  4.  Ion  density  contours  from  various  simulations  to  test  mass  scaling  of  the  wave¬ 
length. 

Figure  5.  Ion  density  contours  from  various  simulations  to  test  Qe/ue  scaling  of  the 
wavelength. 

Figure  6.  Ion  density  contours  from  various  simulations  to  test  scaling  of  the  wavelength, 
varying  mt/me  and  ne/we  to  keep  Rb  constant. 

Figure  7.  Ion  density  contours  from  various  simulations  to  test  scaling  of  the  wavelength, 
varying  VD/c  and  Ue/ue  to  keep  Rb  constant. 

Figure  8.  Comparison  of  simulations  with  different  initial  velocity  distributions:  [top 
Vr  =  ^(constant),  [bottom]  Vr  =  y/2VDr /r0. 

Figure9.  Comparison  of  simulations  with  different  angular  spread  (±0O):  0o  =  180°,  90°,  45 
with  Vr  =  VDr/r0  and  0o  —  90°  with  VT  —  VD. 

Figure  10.  Simulations  at  various  times  (in  reverse  order)  showing  Y-shaped  structure  at 
top  is  not  bifurcating. 

Figure  11.  Simulation  at  various  times  showing  arms  bend  in  the  electron  sense  of  rotation. 

Figure  12.  Simulation  at  various  times  showing  persistence  and  propagation  of  coherent 
density  clumps. 

Figure  13.  Results  of  simulations  showing  the  effect  of  varying  ceil  size:  ion  density 


contours  at  various  times  for  (A)  Ai  =  0.125c/u>e,  (B)  0.0625,  and  (C)  0.03125.  (The 
standard  case  has  Ax  =  0.2). 

Figure  14.  Results  of  linear  theory  for  the  simulation  parameters  of  the  standard  case 
showing  the  growth  rate  as  a  function  of  k  for  the  LHDI  and  the  first  two  harmonics 


of  the  ECDI. 
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